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We present a novel method for the search of high energy extraterrestrial neutrinos in extended regions. The method is based on 
the study of the spatial correlations between the events recorded by neutrino telescopes. Extended regions radiating neutrinos may 
exist in the Galaxy due to the hierarchical clustering of massive stars, the progenitors of all the Galactic accelerators known so far. 
The neutrino emission associated to such extended regions might be faint and complex due to both the escape of cosmic rays and 
the intricate distribution of gas in the environment of the accelerators. We have simulated extended neutrino emission over an area 
of 10° x 10°, where the intensity fluctuations across the region are modelled as a Gaussian random field with a given correlation 
structure. We tested our proposed method over realizations of this intensity field plus a uniform random field representative of 
the spatial distribution of the atmospheric neutrino background. Our results indicate that the method proposed here can detect 
significant event patterns that would be missed by standard search methods, mostly focused in the detection of individual hot spots. 
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i 1. Introduction 

o : 

+3 • The search for the sources of cosmic rays is one of the most 
standing issues in high-energy astrophysics. Major advances in 
the field during the last years have been realized by the consid- 
erable number of instruments dedicated to Galactic cosmic ray 
' ^ studies HI IH SSI 0,12]. Th e on-goingmeasurements of the cos- 
jnic radiation from the Galaxy O, |8[ |90 are helping to construct 
fvQ a comprehensive picture about the sources of cosmic rays, the 
mechanisms of propagation, and the interaction of cosmic rays 
with the gas and radiation fields of the Galaxy 11 1 CHI . Currently, 
t— ( gamma-ray astronomy has a leading role in the exploration of 
the sky at the very high energies. Measurements of the spatial 
and spectral distributions of the diffuse Galactic gamma-rays 
provide information about the propagation and interaction of 



!>■ Galactic cosmic rays 111 ll Il2l Il3ll . whereas localized excesses 
k>| of gamma-ray emission with respect to models of cosmic ray 
j_j propagation can be interpreted as the result of cosmic ray inter- 
C3 actions in regions of enhanced matter density (which are below 
the spatial resolution of current cosmic ray propagation mod- 
els), or due to the presence of Galactic accelerators injecting 
high energy particles which interact close to their sources of 
origin Il4 Il5l llo, Il7ll . The GeV gamma-ray diffuse emis- 



sion is being measured by the Fermi satellite [7], showing that 
gamma-rays are produced throughout the Galactic disk by the 
interactions of cosmic rays with the gas and radiation from the 
interstellar medium, after their diffusion in the Galactic mag- 
netic fields. The Galactic TeV diffuse emission is within reach 



of ground Cherenkov telescopes like H.E.S.S. Ill8ll and Milagro 
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fl^l . instruments which also have provided the first evidence 
of particle acceleration up to TeV energies from a number of 
individual sources ll20ll . 

Further progress is expected in the upcoming future with 
neutrino detectors such as IceCube lElll . In the energy range 
300 GeV < E < 10 PeV, the IceCube experiment at the South 
Pole uses the Earth as a filter to observe neutrinos from the 
whole northern sky. Its large field of view and dynamic range 
offer a unique opportunity to unveil the sites of both produc- 
tion and interaction of the high energy cosmic rays from the 
Galaxy. As in the case of hadronic gamma-rays, neutrino emis- 
sion implies the efficient acceleration of cosmic rays in Galac- 
tic sources, as well as their interaction with matter and radia- 
tion. Surveys of atomic and molecular gas lEiiEill . combined 
with infrared and optical data, have provided large-scale maps 
of the Milky Way's stellar and gas distribution with unprece- 
dented spatial resolution J24ll . despite the uncertainties in the 
distance estimations. From these observations we know that the 
progenitors of all the Galactic accelerators known so far, mas- 
sive stars, are not located at random within the Galactic disk, 
but rather showing a hierarchical structure. 

All present-day star formation appears to take place in gi- 
ant molecular clouds. Massive open star clusters do not form 
in isolation, but tend to be clustered themselves in the so-called 
clusters complexes [25], associated with giant molecular clouds 
following the spiral arms of the Milky Way disk 026il27ll28ll . De 
la Fuente Marcos & De la Fuente Marcos ll26ll found statistical 
evidence of the existence of at least five dynamical families of 
young and massive open clusters in the solar neighborhood (dis- 
tance < 2.5 kpc), which may extend up to 10 degrees in Galac- 
tic longitude. The relevant scenario for high energy astronomy 
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is one in which we can find groups of stellar open clusters in 
an evolutionary stage in which the combined winds of mas- 
sive stars and supernova explosions dominate the energetics of 
the region. This establishes the existence of potential Galactic 
cosmic ray factories in regions of massive star formation [29] 
which may well extend several degrees in the sky, possibly af- 
fecting the spatial distribution of neutrino events recorded with 
neutrino telescopes. 

In this paper we focus on the potential of neutrino tele- 
scopes such as IceCube to detect a significant neutrino event 
pattern associated with these active regions of the Galaxy. This 
could point to the location and distribution of Galactic cosmic 
ray and neutrino sources during the first years of the IceCube 
experiment, when sensitivity might be still too limited for the 
detection of individual steady point-sources. 

A possible neutrino emission in these regions may show an 
intricate, and likely faint, intensity pattern, due to the complex 
structure of the matter in regions harboring massive stars. On 
one hand we have the clustering of potential accelerators, where 
high energy cosmic rays are originated, and on the other hand 
the available target for cosmic ray interactions, consisting of the 
gas left from the parent molecular cloud, which has not been 
used for star formation nor evacuated from the region by the 
radiation and winds from the stars and supernova explosions. 

Both observations and theory about the spatial structure of 
the molecular and atomic gas in a large range of environments 
have shown that it presents some degree of spatial correlation, 
usually characterized in terms of the power spectra ll30l 31 1. 



The causes of the observed spatial correlations are manifold, 
but large scale turbulence and some local effects like the radi- 
ation pressure and winds from the massive stars and supernova 
explosions are the most important S H Si. These ef- 
fects produce spatial fluctuations in the gas density, which in 
turn, may contribute to the fluctuations of a possible neutrino 
signal inside an area populated with cosmic ray sources. 

We offer here an alternative to the search of these poten- 
tial Galactic cosmic ray factories, in a situation in which lit- 
tle knowledge about the high energy processes taking place is 
available. We introduce here the use of analysis of the spatial 
correlations between neutrino events for the discovery of Galac- 
tic cosmic ray sites. We will study the advantages and poten- 
tial of correlation analysis with respect to standard searches for 
neutrino sources. The method presented here and the standard 
search methods are tested on simulated event patterns, in which 
the physical origin of the possible neutrino signal is translated 
into statistical properties of the neutrino intensity field. 

Modelling the exact characteristics of a possible neutrino 
intensity field in a particular region of the Galaxy is out of the 
scope of this paper, and the examples considered aim to illus- 
trate the performance of search methods. Within the region un- 
der study, we assume that the fluctuations in the intensity of the 
neutrino signal vary spatially according to a Gaussian random 
field with a given correlation structure. The last source of neu- 
trinos we consider is the isotropic foreground of atmospheric 
neutrinos, the statistics of which is consistent with a uniform 
Poisson distribution. 

The paper is organized as follows. In section 2 we describe 



the simulations of the event pattern. In section 3 we present 
our proposed method for the analysis of extended regions. Sec- 
tion 4 shows the performance of both our correlation method 
and standard search methods from the results obtained with the 
simulated event patterns. Section 4 has the discussion. 

2. Simulations of the event pattern 

An event pattern consists of a set of event locations r = 
(ri, rj, rpj) within a region (R), at which N events have been 
recorded. In the simplest case, the data set comprises only the 
event locations. However, in some cases we may have addi- 
tional information related to the events which might have a 
bearing with the nature of the analysis. For instance, the ar- 
rival time of the events, for analysis of variable sources ll36ll . or 
the energy of the events for a higher signal to noise discrimina- 
tion power 13711 . These cases correspond to what is known as 
marked point pattern. Although this information can be incor- 
porated in any analysis, in this paper we have not considered 
marks in the event pattern, and we focused only on the spatial 
distribution of events. 

2.1. Simulations of the astrophysical neutrino signal 

A spatial event pattern can be thought of as the realization of 
a spatial process. In many astrophysical scenarios, such physi- 
cal processes depend on random components, and they are usu- 
ally modelled in terms of random fields ll39ll . 

We adopted a Gaussian random field with a given correla- 
tion structure for the intensity fluctuations of a possible neutrino 
signal inside a region that is much larger than the angular reso- 
lution of the experiment. 

If the field is Gaussian, we can determine it completely 
given its mean and correlation function (or power spectrum), 
which facilitates the simulations. The starting point of the sim- 
ulations of Gaussian random fields is the convolution of the 
power spectrum of the field with a white noise (i.e., correlation 
free) random signal ll40ll : 

/(r) = j d 3 ke ik r yfP(k)W(k) (1) 

which gives the intensity of the field at each point r from the 
power spectrum P(k), related to the spatial correlation function 
C(r) by the Fourier transform: P(k) = F~ l C(r). W(k) is the 
white noise, which gives the random amplitude and phase to 
the field. 

If the angular extent of R is large, a spherical harmonic de- 
composition of the field would be needed. In this case, the field 
would have the spectral representation: 



Kr) = ^ ^ a Um Y Um {r) 

1=0 m=-l 



(2) 



where Y[ m is a basis of spherical harmonics, and the coeffi- 
cients ai m are complex uncorrected random variables that de- 
pend only on the covariance or power spectrum of the field Pi. 

In our examples, R is small enough as to assume Euclidean 
space. We generated the Gaussian random fields through Fourier 
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phase randomization by discretizing equation Q] with a Carte- 
sian mesh defined over R. For simplicity, we assumed that R is 
a square region of length L = 10°, and used a mesh of 512 x 
512 nodes, giving a 0.02° sampling interval of R, much smaller 
than the angular resolution of neutrino telescopes. 

The sampling of R is done in Fourier space directly assum- 
ing a power spectrum for the intensity field. The Fourier do- 
main is commonly used in the simulations of random fields be- 
cause the random variables at different points are statistically 
independent, whereas in the spatial domain they have long- 
range correlations that are difficult to sample. 

For simplicity, we assumed an isotropic, homogeneous ran- 
dom field; that is, the correlation between two points r\ and 
depends on the length r = |t| of the vector r = r\ — r%, but 
not on its direction. We adopted a power-law power spectrum 
PQk\) — {a + l^l 2 ) - ' 3 , which has an associated covariance func- 
tion described also by a power-law plus an exponential decay, 
with a,p as the parameters that regulate the decay rate of the 
long range correlations with distance: C(t) oc -fl~ l e~ aT . 

Figure Q] shows a realization of an isotropic, homogeneous 
Gaussian random field with the power spectrum described above 
with a — l,yS = 2. In figure[T]the color scale has arbitrary units, 
but the obtained spatial distribution was used for the represen- 
tation of the astrophysical neutrino signal, where the intensity 
will be normalized to a certain number of signal events in the 
final data sample of neutrino telescopes. The test region over 
which we inject the signal is defined in the range in equatorial 
coordinates 30° <dec< 40°, 300° <R.A.< 310°. 

The instrumentation introduces an additional correlation in 
the intensity pattern due to the errors in the reconstruction of 
the event directions. Therefore, for the analysis of the event 
pattern in section we convolved the signal field of figure Q] with 
the Point Spread Function (PSF) typical of neutrino telescopes. 
The effect of the PSF is to smear the signal over some area, ac- 
cording to its shape. We assume the PSF has a Gaussian profile 
with a = 1, which corresponds to a median of 0.7°, similar to 
what is expected for the IceCube neutrino detector 114 ill . 

2.2. Simulation of the background of atmospheric neutrinos 

A simulated sample consisting only of atmospheric neutrino 
events was constructed injecting uniform random events on a 
sphere. In reality, atmospheric neutrinos reach the detector with 
a zenith dependent component ll42ll . but within a declination 
band of ~ 10° scales like in the examples we consider here, the 
spatial distribution of atmospheric neutrino events in the final 
sample is uniform random, where the intensity at each point 
follows a Poisson process. 

A total of 50000 background events were injected randomly 
through the northern sky, a number which could be represen- 
tative of the total atmospheric neutrino events with E > 500 
GeV expected for the IceCube detector in its final configuration. 
This produces and expectation of 207 background events in the 
10° x 10° region within the declination band 30° < 6 < 40° that 
we have considered in our examples. 
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Figure 1: Realization of a Gaussian random field with power spectrum P(\k\) = 
(1 + |ft| 2 )~ 2 over the 10° X 10° region used in our examples. In this specific re- 
alization units are arbitrary, the intensity of the field is normalized with respect 
to a reference value. In the following examples the signal will be normalized to 
a given number of events. 



3. Analysis of the spatial correlations: the Multi Point Source 
method (MPS) 



Methods to account for the spatial correlations in the data 
are extensively used in astronomy [4311. starting from the work 
of Peebles and co-workers to study the large-scale matter distri- 
bution in the Universe through the two-point correlation func- 
tion of galaxies 114411 - In high energy astrophysics, the searches 
for correlations in the data have been applied mainly to detect 
anisotropics in the spatial distribution of cosmic -rays 1I45I l46ll 
and neutrinos l47ll . 

The most commonly used estimators for the two-point cor- 
relation function in astronomy are based on pair counting (see 
Kerscher et al. (2000) 14811 for a comparison of different estima- 
tors). Following the notation of Szapudi & Szalay (1998) [49], 
the number of event pairs within a distance r is defined as 



P D D(r) = ^^> r (x,y) (x*y) 

xeD yeD 

PRR(r) = £J> r (x,y) (x*y) 

xeR yeR 



(3) 



(4) 



where the summations run over event coordinates in the data 
sample D, and in a sample R of randomly distributed events. 
^VAx, y) = 1 only if a certain pair selection criteria are satisfied 
in the real (random) sample, and is equal to otherwise. After 
the normalization of the number of pair counts: 



DD(r) = ^P DD (r)N(N-l) 
RR(r) = J] PRR(r)N R (N R - 1) 



(5) 
(6) 



with N and N« being the total number of events in the real and 
random sample, the so-call natural estimator of the two-point 
correlation function, rj n , is: 



DD 



(7) 
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which represents the excess probability, with respect to a ran- 
dom distribution, of finding an event at a distance r of another 
event. 

In this section, we introduce this formalism for the detection 
of high energy neutrinos from an extended region. The classical 
two-point correlation function method, briefly described above, 
was adapted to our specific case and optimized for discovery. 
The aim of our method is to determine whether the observed 
spatial properties of neutrino events inside the area under inves- 
tigation are compatible with those expected from a distribution 
of background-only events (atmospheric neutrino events), or if 
an extraterrestrial component has to be invoked in order to ex- 
plain the observations, and at which level the background-only 
hypothesis is rejected. 

Due to the differences our approach poses to the natural def- 
inition of the two-point correlation function, we defined our test 
statistics as a clustering function, and in what follows we will 
refer to our proposed method as the Multi Point Source (MPS). 
Here we describe how the MPS method is defined. 

Let R be the region under study, in which a total number 
of Nmside events have been registered in the final data sample 
at the locations (r\,ri, ...rN„„ ide ), defined by the coordinates r, = 
(0i, (pi) on the sky. MPS makes a two-point sampling of R with 
circular bins of variable area A = n® 2 centered at the locations 
of each of the events inside the region, r,. In the sampling of 
the region, the angular distance, 0,j, from an event i located 
inside R to an event j (at any location) is measured for each of 
the N i nS ide events present in R (z = 1, ...,N inside ; j = 1, ...,N ma i). 
The number of pairs ij is measured as a function of the angular 
separation and the histogram of event pairs as a function of the 
angular distance © is constructed from these measurements. 

The analysis then makes use of a scale dependent cumula- 
tive clustering function, O(0), defined as the excess, with re- 
spect to the null hypothesis, in the number of event pairs within 
a certain angular distance 0: 



<D(0) 



DD(®)d® 



(8) 



f®RR(®)d® 

where DD(@) = Y DD ih RR{®) = V RR U , and the sum runs 

ij U 

over all non-repeated pairs in the real data sample and in the 
background case, respectively. In our case, DDjj(RRjj) = 1 only 
if either the event i or the event j, or both, are within the region 
under study, and it is equal to zero otherwise. That is, MPS 
does not correlate events with specific locations, but rather it 
considers each event inside the region as a point source in order 
to determine its degree of correlation with the rest of the events. 
This feature is precisely what names the method as "Multi Point 
Source". With this definition we measure both the intensity of 
the process that generated the observed neutrino event pattern 
as well as its correlation structure. 

The ultimate result of the analysis in neutrino telescopes is 
the p-value of the observation, that is, the probability that the 
observed event pattern is just a realization of the atmospheric 
neutrino background field. A discovery is claimed when this 
probability is below 2.8 x 10~ 7 , which corresponds to a 5<x 



detection. To obtain the p-value of an observation, we must 
therefore know the distribution of our test statistics under the 
background-only hypothesis. In the MPS method, the test statis- 
tics is a function of the angular scale 0. For each exam- 
ined, we must construct the probability distribution of <D(0). 
The number of background events in a given area A = n® 2 
is produced by a Poisson process; as a consequence, the num- 
ber of pairs of background events within a distance follows a 
Gamma distribution. However, with the definition of event pairs 
in MPS, the parameters of the Gamma distribution are difficult 
to estimate from observables, and we obtained the distributions 
with Monte Carlo simulations. 

4. Analysis of the event pattern 

In this section, we present our results about the ability of 
both a stardard search method and the MPS method to find a 
significant neutrino events pattern over extended regions. The 
tests were performed over simulated skymaps, in which a region 
R of 10° x 10° exhibits neutrino emission with some correlation 
structure (section l2~TT i. 

Standard search methods deal with the identification of hot 
spots, defined as a spatial concentration of events. The search 
for hot spots is done by mapping the sky using different tech- 
niques; the ones commonly used in neutrino astronomy are 
based on either "classical" binned methods, or likelihood-based 
methods with different degrees of sophistication in the mod- 
elling of the data 1 37, 3^1. Regardless of the method used in the 
analysis, the searches have been usually optimized for the dis- 
covery of point-like emission, and hence correlations of events 
beyond the single point-source scale do not enter in the analysis. 
The exploration of the sky, seeking for sources of high energy 
neutrinos, is done by superimposing a grid over the event dis- 
tribution. Classical binned methods count the number of events 
within a circular fixed-area bin centered at each of the nodes 
in the grid, resulting in estimates of the intensity of the pro- 
cesses that generate the observed event pattern at each sample 
point. The grid step is much smaller than the angular resolution 
achieved in the analysis, and the bin size is usually optimized 
to have the best signal to noise ratio for point-sources. Likeli- 
hood methods use the same type of scan analysis, but they give 
a probabilistic assignment of an event to a component (signal 
or background) according to its distance to the point at which 
the intensity is being estimated. At each grid point, the likeli- 
hood of a mixture model of signal plus background is compared 
to a pure background hypothesis, and the region under study is 
imaged as a probability density map. 

The performance of both classical and likelihood methods 
are similar in terms of discovery, as long as they use the same 
event information 13711 . Given that the event density maps ob- 
tained from a binned analysis offer a clear visual inspection of 
the spatial event pattern, as well as a more straightforward com- 
parison with the results from MPS, we used a classical binned 
method as the standard scan analysis of our simulated skymaps. 

The scan of the test region of 10° X 10° was done in steps 
of 0.25° x 0.25°. At each point, the local event density was 
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measured and compared to the one expected under the back- 
ground hypothesis, where by local we mean within a fixed-area 
bin around the sample points. 

Both the MPS and the binned scan of R were applied to 
background+signal samples and the corresponding significances 
were calculated with the p-value of the outcomes of each method. 
The probability distributions were computed after the analy- 
sis of n random background-only samples with both the MPS 
method and with a binned scan of R. Each of these simulated 

(n) 

data sets yields n event patterns with a different number N. ., 

J 1 inside 

of background events inside R. The average clustering function 
in the background case and its dispersion at different angular 
scales is obtained after applying MPS over n — 10 4 background- 
only data sets. The obtained probability distributions were fitted 
to a Gamma distribution. The binned scan was applied to the 
same background-only datasets and the distribution of the num- 
ber of events within the search bin for each of the nodes in the 
grid was fitted to a Poisson distribution. 

In the next sections, we present our results under two differ- 
ent scenarios: a situation of high signal to noise ratio (S/N), and 
a situation in which the S/N is very low. The case of a high S/N 
underlines the differences in the information extracted from the 
same event pattern with the two different methods used in this 
paper. The study of a low-signal case illustrates how the MPS 
and a standard scan of R deal with the background dominated 
samples characteristic of neutrino telescopes. 

4.1. High S/N case 

In this example a total of 100 signal events were injected 
randomly following the distribution of figure [T]convolved with 
a Gaussian with cr — 1, representative of the PSF of current 
neutrino telescopes. On top of the signal, background events 
were randomly generated as explained in section l2~2l 

Fig. [2] shows a realization of the signal+background field 
imaged as a significance event density map from the results 
of the binned scan of R with circular bins of radius 0.5°, 1°, 
1.5°, 2°. In this simulated data set, a total of Ni„ S id e = 294 
signal+background events have fallen within R, yielding a 6<x 
excess of events inside the entire R, clearly representing a high 
S/N case. 

The appearance of the resulting images is a combination of 
diverse effects: the spatial structure of the signal intensity field, 
the smearing of this signal because of the instrumental error in 
the event locations, and the size of the search bin area. The sig- 
nificances obtained then reflect the balance between the clus- 
tering of the signal at the scales of the search bin radius, and 
the scales at which random accumulation of background events 
is more probable. A discovery (i.e. a p-value < 2.8 x 10~ 7 ) is 
achieved when sampling the region with bins of 1 .5° radius (fig- 
ure |2(c)| i, but the significance drops considerably when observ- 
ing at different angular scales, and any information regarding 
the spatial structure of the signal within the region under exam- 
ination is lost. However, this information is recovered when we 
study the spatial correlations between events with MPS. The re- 
sults of the analysis of the region with MPS are shown in figure 
[3] Figure |3(a)| shows the clustering function of the events up 
to 5° scales, and figure [3(b)| shows the significances obtained at 



each angular scale. The measured spatial correlations between 
events below 5° follow a power-law with exponential decay, 
consistent with the signal field that was used in the simulations, 
and the discovery of the simulated astrophysical neutrino sig- 
nal goes up to angular scales of ~ 2.5°. Therefore, MPS results 
have the potential to discover a more complicated structure of a 
possible astrophysical neutrino signal than what was observed 
with the standard scan of the region. 

4.2. Low S/N cases 

In this example a total of 40 signal events were injected 
in the region. With a background expectation of 207 events 
within R, this situation corresponds to a background dominated 
sample. This case exemplifies a more realistic situation, since a 
total of 40 signal events present at the final sample of neutrino 
telescopes from a region of 10° x 10° constitutes a more realistic 
situation given the current effective areas of neutrino detectors 
and sources with moderate neutrino fluxes ll4ll . 

We illustrate here three simulated examples of such situa- 
tion. Fig. [4] shows the result of a 1° binned scan for three real- 
izations of the background field plus 40 signal events following 
the spatial distribution of figure Q] The three cases are compat- 
ible with background fluctuations, reaching a 3<x level at most. 
Fig. [5] shows the same scan, but in which signal events have 
been removed. The comparison of figs. |4] and [5] manifests the 
random appearance of events even in the signal+background 
maps, where the most significant spots occurred always at ran- 
dom locations, as expected in background-only skymaps, but 
with a slightly higher significance due to the additional pres- 
ence of signal events in the region. If we combine the three 
simulated samples (figures |6(a)||6(b)| i, corresponding to a larger 
period of data taking in the real situation, the S/N is still too low 
to produce a significant detection from the binned scan, and the 
minimum p-value achieved only reaches 3cr at 0.5° scales (fig- 
ure |6(a)|. If, in addition, we correct from the trials associated 
to the scan of an extended region, this significance is further 
reduced. However, although the significance of these individ- 
ual spots is not enough for a detection, their distribution is not 
common of a random distribution of background events. This is 
illustrated when we use MPS to quantify the departures of the 
whole event pattern within R from a distribution of background 
events. With MPS, the minimum p-value of the observed event 
pattern reaches 4.3cr observation (figure |6(c)| i. The trial fac- 
tors from scanning several clustering scales with MPS are very 
small, and the MPS post-trial significance remains practically 
unchanged, from the pre-trial p-value of 10 -5 5 to a post-trial 
p-value of 10~ 5 ' 3 , yielding also a 4.3<x final result. 

Therefore, in the typical situations of neutrino telescopes of 
low S/N, the only hint of a possible neutrino signal would be the 
accumulation of relatively weak spots within the region which 
is believed to be a factory of cosmic rays. As long as the S/N 
is low, these most significant spots do not necessarily have to 
appear at the same location in the analysis of the same region 
in different epochs, since they have an important background 
component, which fluctuates randomly within the region. The 
signal would be missed by standard scan analysis, but could 
be detected by correlation analysis like MPS. Eventually, after 
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(a) Significance map of event densities at 0.5° scales 
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(b) Significance map of event densities at 1° scales 
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(c) Significance map of event densities at 1.5° scales 
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Figure 2: Results from the binned scan. The color scale is in units of standard deviation. 
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Figure 3: Results from the MPS 



some years of integration time, the S/N would be high enough 
to overcome the background fluctuations, providing an image 
of the neutrino emission in the region. But before this situation 
is achieved, we expect to find a distribution of neutrino events 
which departures significantly from the average distribution of 



background events, provided there is some clustering of neu- 
trino sources within the region under examination. 
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(a) Significance map of event densities at 1° 
scales in sample 1 
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(b) Significance map of event densities at 1° 
scales in sample 2 
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(c) Significance map of event densities at 1 ° in 
sample 3 



Figure 4: Results from the binned scan over signal+background samples 
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(a) Significance map of event densities at 1° 
scales in sample 1 
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(b) Significance map of event densities at l c 
scales in sample 2 



(c) Significance map of event densities at 1° in 
sample 3 



Figure 5: Results from the binned scan over background-only samples 
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(a) Significance map of event densities at 0.5° 
scales in the combined sample 
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(b) Significance map of event densities at 1° 
scales in the combined sample 



(c) Significance of the clustering function at each 
of the angular scales tested in the combined sam- 
ple. 



Figure 6: Results from the binned scan and from MPS over the combined sample 



5. Discussion 

The discovery of high energy neutrino emission from a par- 
ticular region of the Galaxy would suggest a local injection of 
high energy protons or heavier nuclei; as well as their inter- 
action not very far away from their sources of origin, before 
they lose their energy by diffusion in Galactic magnetic fields 
|5Cj|5l|]. Whereas the discovery of individual point sources of 
high energy neutrinos seems challenging because of the limited 
sensitivity of neutrino detectors, the agglomeration of neutrino 
sources inside the same region of the sky could be detected as a 
significant departure from the distribution of atmospheric neu- 
trino events. Such clustering of neutrino sources would be a 
natural consequence of the distribution of stellar matter in the 



Galaxy. Massive and young stellar populations, the birth sites 
of Galactic core-collapse supernova, are associated to the spiral 
arms of the Milky Way, following a hierarchical scheme where 
stars are grouped in clusters, and the stellar clusters themselves 
form part of open cluster complexes. 

Under the assumption that the explosion of a supernova 
is the most plausible energy source for cosmic ray accelera- 
tion in the Galaxy, we propose in this paper to search for neu- 
trino emission in the regions where massive stars die, that is, 
where the cosmic ray acceleration process takes place in the 
supernova remnant shocks that form around the explosion site. 
Given the distribution of massive star clusters in the Galaxy, we 
can find such potential cosmic ray factories at angular scales 
which could even reach around ~ 10° or more, depending if we 
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also want to consider clustering of potential accelerators due 
to projection effects. This type of search therefore implies the 
analysis of extended regions, much larger than the detector's 
angular resolution. We emphasize the usefulness of a search 
method which accommodates different degrees of complexity 
in the spatial distribution of astrophysical neutrino events, and 
try to not restrict ourselves only to the case of point-like emis- 
sion. In this paper, we propose a method to search for neutrino 
emission from extended regions which takes into account all 
these considerations, and that we conveniently named the Multi 
Point Source (MPS). Comparing to a single source search, the 
MPS would give a better result only if the emission deviates sig- 
nificantly from spherical symmetry, whereas it would be com- 
parable to a standard binned search in the opposite case. 

In situations in which the search for extraterrestrial neutri- 
nos goes beyond the single and spherically symmetric source 
approach, we have shown that correlation analysis like MPS are 
able to discover a significant neutrino event pattern inside the 
region under investigation with a higher significance than stan- 
dard searches, which scan the region looking for single sources. 

In simulated examples of a high signal to noise ratio, where 
also a standard scan of the region results in a discovery, the 
MPS is proven to be worth performing in order to extract infor- 
mation of the neutrino event pattern at different angular scales, 
where no significant information is extracted from the scan. 
Even more interesting, given the current status of neutrino as- 
tronomy, is the case of a very low S/N, where we have shown 
that even when the standard scan of a particular patch of the sky 
do not yield any hint of an astrophysical signal, MPS indicates 
the existence of astrophysical neutrino signal inside the region 
under study. 

The conclusions of this paper about the potential of the 
search method we propose to discover neutrino emission with 
respect to standard search strategies are independent on the spa- 
tial correlation function that we adopted for the simulations. 
The examples used here are hypothetical, and not attempt to de- 
scribe the neutrino emission from any particular region of the 
Galaxy, but to illustrate the usefulness of our method for de- 
tecting a significant event pattern even in strongly background 
dominated samples. 

The authors are thankful to Aday Robaina for useful discus- 
sions about the two-point correlation function. Yolanda Sestayo 
acknowledges the financial support from J.M. Paredes through 
ICPvEA Academia. 
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